Characterisation of spatial network-like patterns from junctions' geometry. 
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We propose a new method for quantitative characterization of spatial network-like patterns with 
loops, such as surface fracture patterns, leaf vein networks and patterns of urban streets. Such 
patterns are not well characterized by purely topological estimators: also patterns that both look 
different and result from different morphogenetic processes can have similar topology. A local 
geometric cue -the angles formed by the different branches at junctions- can complement topological 
information and allow to quantify the large scale spatial coherence of the pattern. For patterns that 
grow over time, such as fracture lines on the surface of ceramics, the rank assigned by our method to 
each individual segment of the pattern approximates the order of appearance of that segment. We 
apply the method to various network-like patterns and we find a continuous but sharp dichotomy 
between two classes of spatial networks: hierarchical and homogeneous. The first class results from 
a sequential growth process and presents large scale organization, the latter presents local, but not 
global organization. 



PACS numbers: 89.75.Fb 89.75.Hc 47.54.-r 
I. INTRODUCTION 

A central issue in complex systems research is to un- 
derstand the formation and the properties of spatio- 
temporal patterns found in physics and biology. To this 
end, an essential step is the definition of appropriate mea- 
sures of their form: quantitative estimators that allow 
to compare different patterns and to objectively validate 
models. This paper focuses on the family of spatial pat- 
terns that are "network-like" . 

Network-like patterns are common in natural and ar- 
tificial systems: they are found in leaf veins, fractures 
on the surface of materials, patterns of urban streets and 
animal trails / galleries, river networks, blood vessels and 
circulatory systems. The factors underlying the forma- 
tion of these patterns are different from system to sys- 
tem: surface cracks result from the shrinkage and stress 
of materials; urban streets are generated by human ac- 
tivity and so on. In spite of intrinsic differences, a few 
simple morphogenetic events describe the formation of all 
these patterns: nucleation of new network components, 
elongation of existing segments, branching and intersec- 
tion. The final topology of the pattern is completely 
determined by the sequence of such growth events, plus 
a pruning event: the cut or removal of already formed 
segments. 

Many studies have used topological estimators to de- 
scribe the form of network-like patterns and better un- 
derstand their morphogenesis and functional properties. 
One of the first progresses in this direction goes back 
to the Horton-Strahler coefficients 0, HI, a simple but 
powerful tool to describe quantitatively the form of hy- 
drogeologic networks. The coefficients, which work by 
assigning rank numbers to every branch of a tree-like net- 



work, were successfully applied to the study of different 
systems, from river networks to leaf patterns [f|, 

and even to ant trail 0] and termite gallery [8[ patterns. 
Unfortunately, Horton-Strahler's coefficients cannot be 
computed on networks with cycles. This has a profound 
impact on the applicability of the method, as many real- 
world network-like patterns have cycles. 

General graphs are characterized by a whole set of mea- 
sures of the local and global organization: node degree, 
assortativity between nodes, clustering coefficient, fre- 
quency of specific subgraphs, presence and number of 
cycles, diameter, path length etc. Several studies 

have computed graph measures and used them to quan- 
tify the local and large scale organization of real world 
spatial patterns, such as urban street patterns fl2i-[20lj . 
systems of animal trails [2lJ and galleries 123426 1. net- 
works of channels in trabecular bones [U l28| |. networks 
of fungal hyphae [29j . 

Unfortunately, purely topological estimators, however 
useful, do not fully account for the organization of 2D 
network-like patterns which undergo specific constraints. 
For instance, in planar graphs, the average node degree 
cannot be higher than six [llj . In addition, for a number 
of network-like patterns, the degree distribution is even 
more regular than what is imposed by planarity: the 
large majority of nodes is found to have degree equal to 
three in very different systems, such as two-dimensional 
foams [HI, honeycomb patterns and biological epithe- 
lia [33[, leaf vein networks [Hj], fracture patterns [35j . 
In self-organized, "bottom-up" towns [lH, [l5[ the ma- 
jority of junctions between streets also have degree three 
(though the same does not hold for planned towns, where 
the majority of crossroads have degree four [l5[ or even 
higher [361 ] ) . The homogeneity of degree is reflected by an 
homogeneity in the length of network cycles that are com- 
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FIG. 1: (Color online) Drawing of the evolution of a fracture 
pattern in a solution of latex particles [37]]. In each image, the 
yellow (light gray) lines indicate the newly appeared fractures. 
When new lines of fracture appear their growth is affected by 
the already existing fractures (for instance, newer fractures do 
not cross older ones). Conversely, the position of older lines 
does not change when new ones appear. 



posed by six edges on average in all the aforementioned 
systems. Similarly, network distances are not very infor- 
mative as they basically scale with euclidian distances in 
all these systems. 

Here we want to use the spatial information to provide 
a deeper understanding of the pattern. We focus in par- 
ticular on the information carried by the angles formed 
at junctions. This choice is motivated by two comple- 
mentary arguments, one about the physics of the growth 
of the network-like pattern, which can influence directly 
the branching angles, and the other about dynamic prop- 
erties (e.g. navigation) that may be associated with the 
resulting network. 

Throughout the text, as we consider spatial networks, 
we will make extensive use of the words edge, arc and seg- 
ment with the following definitions: an edge is a curve 
-generally almost linear- between two consecutive junc- 
tions (classical network edge); an arc is a directed edge 
when an orientation is added to the network; a segment is 
a series of contiguous edges grouped together as described 
later in section Inland used to reveal a structuration at a 
larger scale. 



A. growth 



segments is affected by the older segments; the older seg- 
ments, however, do not undergo further reorganisation 
as a consequence of the appearance of new ones. Under 
these conditions (sequential process, growth with no dele- 
tion of edges and absence of reorganization) the temporal 
hierarchy in the appearance of segments is reflected into 
the spatial hierarchy of lengths and arrangements of the 
final pattern [33] ■ Across a junction the edges belong- 
ing to the older segment are the straight continuation of 
one another; conversely, the edge belonging to the newer 
segment typically forms a large angle with the others. 



But the local angle information is not sufficient in itself 
to attribute the edges to the same or a different segment: 
consider the case of a junction like the one in figure [2} A. 
In the absence of any information, the most natural infer- 
ence is to group edges e± and e2 into the same, older, seg- 
ment. The attribution is more problematic if one already 
knows that e^ is older: then it would be more natural to 
group e3 and e\ into the same old segment, and assume 
e2 appeared later. In summary, the best inference about 
what network edges should be grouped into the same seg- 
ment depends on the continuity between adjacent edges, 
but also on the information already available about the 
order of appearance of other segments. It can be con- 
structed only sequentially. 





e! 



Network-like patterns often are formed as the result 
of a sequential process, where new segments appear at 
different times [35] without undergoing further reorga- 
nization. One example of such patterns is provided by 
the fracture lines on t he g laze of ceramics illustrated in 
figure [U Bohn et al. [35J suggested that urban streets 
patterns fall into this same morphological class. 

The appearance, elongation and termination of new 



FIG. 2: (Color online) A: Small graph, involving a sin- 
gle junction. B: The weighted directed line graph built 
from the graph in A. Weights are as follows: red (dark 
gray) — > weight=0; green (light gray) — > weight=l 
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B. navigation 

When the network-like pattern is also the support for 
a transportation function, as it is for instance the case 
with patterns of urban streets, the angles of edges at 
intersections and junctions play a role for navigation and 
orientation. 

In architecture, this concept has been theorized under 
the name of "space syntax" (38|. With the space syn- 
tax method of analysis, urban street patterns are frag- 
mented into a number of straight segments (ideally the 
minimum number of segments in which the line of view is 
conserved) and one analyses the network obtained repre- 
senting straight segments as nodes, and intersections be- 
tween segments as edges. Then selecting a line (a street) 
as a starting point, one can number each line in the map 
according to how many changes of direction separate it 
from the starting line. This measure is generally referred 
to as "depth" and is a measure of distance: it repre- 
sents the minimum number of changes of direction to go 
from the origin to any other place in the street network. 
This kind of measure has been correlated with different 
aspects of social life: patterns of pedestrian movement, 
traffic, property value and so on ;39]. Nevertheless, the 
process of fragmenting street patterns into segments with 
the same line of sight was not proven to always have 
a computable solution and it might be too sensitive to 
small differences of orientation of the urban grid [401 ] . 

In the physics literature, analogous ideas are followed 
through what is usually called a "dual network" approach 
(though a more appropriate definition would be "line net- 
work" approach) , where named streets [HI, [U |42[ or con- 
tiguous segments [l|| are mapped into network nodes and 
there is an edge whenever two streets intersect or bifur- 
cate. Dual networks also provide a convenient represen- 
tation for measuring the amount of information neces- 
sary to navigate inside the network j4l| , in particular for 
navigation strategies relying more upon the continuity of 
linear elements than on salient points. 

In fact one can imagine to describe a path through the 
network with instructions of the form: "go straight for 
N steps (N junctions or crossroads), or until you find a 
salient point (a traffic light, a square...)", followed by in- 
formation on which new direction to take. In general a 
"simple" path will be one that involves only few devia- 
tions from the current direction, not necessarily a short 
one. However, the result will not be the same depending 
on the direction followed. Let's illustrate this with ref- 
erence to the junction of figure [2]-A. Moving from edge 
e\ to e%, or in the opposite direction from e2 to e\ does 
not involve any change of direction. Going from e\ to 
e3 requires one change of direction (i.e. to abandon the 
straightest path and take an alternative one) , but the op- 
posite is not true since e\ is the more direct continuation 
of movement when coming from e%. 

Both the argument about network growth and the ar- 
gument about navigation illustrate the interest of classi- 
fying the edges of network-like patterns using information 



provided by junctions' angles. Our exaxmples also indi- 
cate that the process is not symmetric with respect to 
the starting point used for classification. 

In the following section [TT] we introduce a new algo- 
rithm that given a local rule (namely an evaluation of 
the angles at each junction) and a set of arbitrarily cho- 
sen root edges, assigns a rank (expressed by an integer 
number) to each edge of the network-like pattern. Con- 
tiguous edges with the same rank can then be grouped 
together into segments. The numbers of segments for 
each given rank and their average properties provide a 
quantitative statistical characterization of the morphol- 
ogy of the pattern. Possible applications to lattices and 
real world spatial networks are explored in the remaining 
sections Hill to IVlTl 



II. ALGORITHM DESCRIPTION 

The classification and numbering of segments involves 
two steps. First, given a root edge, or a set of root edges, 
a rank number is assigned to each edge. Then, the con- 
tiguous edges with the same rank are grouped into seg- 
ments. 



A. computing edge ranks 

Let us denote by G a graph which models a network- 
like structure. Given a root edge e r in G, the rank of an 
edge ei of G expresses the minimum number of direction 
changes needed to reach when coming from e r . We 
can also take into account a set S r of root edges; in this 
case, the rank of e^ is the minimum number of direction 
changes when coming from the closest root edge in S r . 

The propagation of ranks across a junction depends 
on the direction in which the junction is traversed. For 
instance, in the graph of figure [2J-A if we set the rank of 
edge ei to zero (r (ei) = by convention), then the rank 
of e3 is equal to 1 (r(ea) = I), but the same relation does 
not hold in the opposite direction: if we set the rank of 
63 to zero (r(e3) = by convention), then the rank of e\ 
is also equal to zero (r(ei) = 0: no deviation from the 
straightest path when going from e% to ei). 

Generally speaking, each edge e^ can be crossed in two 
directions; by convention, we denote ef (resp. ef ) the 
arc associated to ei when it is crossed towards the right 
(resp. towards the left). To make the computation of all 
the ranks easier, we consider the directed line graph Go 
(fig [2J-B) which models the connections between all the 
arcs. The vertices of Go are the arcs {ef 1 , ef} associated 
to any edge of G, and two vertices e* and e* (the star 
means either R or L) are linked together in Go if the 
destination of e* is also the origin of e*. 

For instance, the crossing of the path e\ — e2 from the 
left is modeled by ef — > ef in Gd, and the crossing in the 
opposite direction by ef — > ef . There is no connection 



between e\ and as the navigation e\ — > is not 
directly possible. 

Weights are assigned to the arcs of Gd as follows. The 
weight w (e*, e*) is equal to when the path formed in 
the initial graph by the edges and is the straightest 
one and 1 otherwise. We can also introduce a threshold 
on the maximum angle, implying that if the straightest 
pairing of edges still involves a too large angle, then it 
also will have weight equal to 1. Formally, w (e*, e*) = 
when (1) the angle formed in the initial graph G by e* and 
e* is the minimum angle formed by e* and its adjacent 
edges, and (2) this angle is smaller than a given threshold 
(45 degrees in our analyses). 

Once defined a root edge e r in the initial graph whose 
rank is equal to by convention, the rank r (e^) of any 
other edge can be easily obtained by a distance com- 
putation in Gd' 

t fe) = 

Mm {d (e* ef ) , d (ef, ef ) , d (e? , ef ) , d (ef ; , ef ) } 

where d is the shortest path length in the directed line 
graph Gd weighted by w. When a set S r of root ver- 
tices is considered, the rank of is the minimum of the 
ranks computed from all the edges of S r with the above 
formula. 



B. grouping edges into segments 

Several adjacent edges with the same rank can be 
grouped together into "segments" . A segment is a series 
of edges all having the same rank and with exactly two 
ends (not including junctions). The segment is initial- 
ized with any edge of G. Its continuation is determined 
by pairing together adjacent edges of the same rank until 
one of the following conditions is encountered: 

1 . There are no more edges with the same rank at one 
extremity. 

2. The segment intersects an edge of lower rank. 

In the special case when three or more edges with the 
same rank are incident to the same junction, the contin- 
uation of the segments is determined by pairing together 
the two edges that provide the straightest continuation 
of one another. Then, we proceed in a similar way pair- 
ing together all the other edges incident to the node. If 
an odd number of edges with the same rank meet at a 
junction, one edge remains excluded from all the pair- 
ings and the segment containing that edge is ended at 
the junction. 

We can give an intuitive justification for the above 
rules. Imagine the case of a network-like pattern re- 
sulting from a sequential growth process, such as for in- 
stance fracture lines on ceramics. Grouping edges into 
segments is equivalent to identify those network edges 
that belong to the same line of fracture. With condi- 
tion [5] we impose that younger lines of fracture (higher 
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FIG. 3: (Color online) Illustration of the procedure for 
grouping together edges into segments. Al: the red (near- 
horizontal) edges have lower rank than the yellow (near- 
vertical) edges. The edges are grouped in A2 into three seg- 
ments: a green (near-horizontal) lower rank segment and two 
(black and cyan) segments with higher rank. Bl: all the edges 
have the same rank. In this case, two straight segments are 
identified in B2. If an odd number of segments with the same 
rank meet at a junction (CI), the one that deviates more from 
the direction of the others forms a distinct segment (C2) 

ranks) never cross already formed fracture lines. A paral- 
lel can also be found with what happens for urban street 
patterns, where small streets usually change their name 
when crossing larger ones. 

Figure |3] shows examples of how segments are grouped 
together. In Al the near-horizontal edges have lower 
rank than the near- vertical edges. In this case, the near- 
vertical edges form two distinct segments (A2), in spite of 
being in direct continuity. Bl presents exactly the same 
pattern, except that now all the edges have the same 
rank. In this case, the edges are grouped together in two 
intersecting segments, one near-horizontal and one near- 
vertical (B2). When there is an odd number of edges 
incident to a junction, as in CI, the couple(s) of edges 
that deviate less from each other's direction are paired 
together in the same segment(s) and the remaining edge 
forms a segment by itself, which is ended at the junction 

To summarize, the straight contiguity of edges is well 
represented by a directed line graph Gd whose arcs are 
weighted with appropriately chosen weights. The arbi- 
trary selection of a set of root edges in G allows to ex- 
press changes of direction in terms of a distance mea- 
sure, allowing to assign each edge of G a ranking number. 
Edges that are more likely to belong to the same spatio- 
temporal event of network formation can be grouped 
together into segments. The C++ computer code im- 
plementing the rank computation and segment mapping 
is available as online supplementary material, alongside 
with sample networks to test it. 



III. LATTICE MODELS 

In this section we explore the distribution of segment 
ranks and edge ranks in a simple lattice model. In par- 
ticular, we introduce a class of lattice models that we 
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FIG. 4: (Color online) Mondrian pattern generated by iter- 
ated domain divisions. Here, divisions are alternating: all the 
cuts made at division t + 1 are orthogonal to the cuts at di- 
vision t. In all figures, the color (shade of gray) of the edges 
reflects their rank; a periodic rainbow colormap is used. For 
this example, all the edges in the bottom of the figure are 
selected as roots. 



call "Mondrian" lattices (figure U]) , intended to mimic 
the growth process and the characteristics of hierarchi- 
cal network-like patterns. In order to build the lattice, 
we start with a single rectangular domain and we iterate 
the following operation: each rectangular cell of the lat- 
tice is divided in two smaller rectangles by introducing a 
new cut parallel to one of its sides (let's say horizontal 
or vertical) at a random position chosen from a normal 
distribution around its centre. The cuts can either be 
"alternating", where horizontal cuts at time t are fol- 
lowed by vertical cuts at time t + 1, and vice- versa, or 
"random" , where the horizontal or vertical direction of 
cut is chosen randomly for each cell and each division. 

Here, and in the rest of the paper, we will focus mainly 
on two kinds of statistics. The first is the histogram of 
segment ranks (the frequency of occurrence of segments 
with a given rank). The second measure relies on edges 
(not segments) and looks at how the ranks of individual 
edges increase with their topological distance from the 
roots. The topological distance between an edge and 
a root e r is the minimum number of junctions on a path 
between e r and in G. Intuitively, the rank measures 
the number of direction changes needed to reach a par- 
ticular edge and the topological distance measures the 
number of vertices crossed. Hence, if we plot the rank 
r(ei) of each network edge against its topological dis- 
tance from the roots td{ei), we obtain a characteristic 
distribution that can be fitted with a linear function of 
the form 



r(e 4 ) = a ■ tdfe) + f3. 



(2) 



The slope of the fitting line is tg(a), and its inverse 
, gives the average length of the straight segments. 



This can also be interpreted as the largest scale at which 
spatial organization is observed. 

Figure [5l top row, displays the histogram of segment 
ranks for Mondrian lattices after 6, 12 and 18 iterations 
of divisions. The orange and cyan histograms are for 
alternating and random divisions, respectively. For ran- 
dom divisions the asymptotic distribution can be approx- 
imated by a Normal distribution: 



where 



N(r) = 

2+ 



A-e" 



(3) 



6 -J (where t is the number of it- 
erations of the division process), A — \ +1 ■ 1-043 and 

(7 = 1.020--^j±-i. (Details on how the parameters A, r and 
a were estimated are provided in the appendix). In prac- 
tice, for a pattern observed at a given time, t is usually 
unknown, and it is easier to approximate it in terms of 
the number of cells C, assuming that the relation C — 2* 
holds . The number of cells for a planar graph in turn 
is obtained from the number of edges and vertices in the 
network through the Euler's formula, V — E + C = 1 
(where V is the number of vertices and E the number 
of edges). In this way we can compute the theoretical 
rank distribution for lattices of given size. Such theoret- 
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FIG. 5: (Color online) Top row: histogram of segment ranks 
for Mondrian lattices of different size (from left to right the 
number of cells (C) is equal to 2 6 , 2 12 , and 2 18 ). Orange 
(dark gray) histogram: alternating divisions. Cyan (light 
gray) histogram: random divisions (in the case of random 
divisions different realisations are possible: the whiskers on 
top of the histogram give the standard error of the mean). 
Bottom row: distribution of the rank vs. the topological dis- 
tance of lattice edges for the same "alternating" Mondrian 
lattices. Red continuous line: linear fit to the data of the 
form r(ej) = a ■ td(ei) + /3, with td(a) the topological dis- 
tance of edge d from the root edges. The slope (tg(a)) is 
0.148 for the lattice with C = 2 6 , 0.018 for the lattice with 
C = 2 12 and 0.002 for the lattice with C = 2 18 . 
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FIG. 7: (Color online) From left to right: same plots as in 
figure [5] corresponding to three lattices with 2 12 cells each, 
but different levels of degradation (0, 70 and 140 iterations of 
degradation respectively). When the local junction geometry 
is destroyed, the histogram of segment ranks (plots on the 
top) shifts from being peaked to uniform. The fitted slope of 
the rank vs. distance plot (tg(a)) (bottom plots) is 0.014 for 
A 0.355 for B and 0.707 for C. Small differences between the 
graphs in A here and the corresponding ones in the middle col- 
umn of figure [5] are due to the additional nodes introduced by 
imposing periodic boundary conditions during network degra- 
dation. 



FIG. 6: (Color online) The same lattice as in figure [4] where 
the position of each node is shifted by associating to all the 
adjacent edges a vector of unit length and centrifugal direction 
and computing the vectorial sum over all the adjacent edges. 
The pattern in B and C correspond to different numbers of 
iterations of degradation. 



ical distributions are reported as dotted lines on top of 
the histograms of figure[5l top row, as well as in different 
figures throughout the paper. 

Figure [5j bottom row, plots the rank r(e<) vs. the 
topological distance tdfe) for each network edge e, of 
alternating Mondrian lattices of different sizes (6, 12 and 
18 iterations of division), together with a linear fit of 
the distribution (continuous red line). We can see that 
when the number of iterations of division increases, the 
overall slope decreases, revealing in fact a curved relation. 
This is not surprising, as topological distances typically 
scale with the square root of the lattice size, while ranks 
increase with the logarithm of the size. 

These simple lattice models may not be representative 
of the configurations of real world patterns with less reg- 
ular edge orientations. In order to assess the robustness 
of rank indices when the edge orientations are altered, we 
proceeded as follows. Starting from a lattice of 2 12 cells, 
obtained through alternating divisions, we iterated the 
following operation: at each time step the node positions 
are shifted by computing the resultant of vectors directed 



along the incident edges and oriented centrifugally. Each 
vector has magnitude arbitrarily fixed to 30 qq , where L 
is the length of the side of the grid. This progressively 
transforms the Mondrian lattice of figure @JA into the 
foam-like patterns of [BJ-B and C. In order to prevent the 
lattice from shrinking during this process, we impose pe- 
riodic boundary conditions by merging together the op- 
posite sides of the lattice. After the degradation phase, 
the boundaries are separated again, and the rank com- 
putation proceeds normally, except that new nodes have 
been introduced along the perimeter of the lattice, dou- 
bling the nodes already existing on the opposite side of 
the lattice. 

Such an evolution from orthogonal to symmetric junc- 
tions is not a purely geometric artifice: many real world 
patterns do actually evolve from an orthogonal to a sym- 
metric configuration, as the result of forces that mini- 
mize the local length of the network edges. For instance, 
in basaltic lava rocks, the junctions, initially appearing 
at right angles, progressively evolve towards a symmet- 
ric configuration with angles of 120 degrees as the joints 
grow inward during solidification of lava [IH In a 
similar way, leaf veins seem to intersect at right angles 
when they first appear, but then to evolve continuously 
towards more symmetric foam-like junctions 



34] 



Figure [BJ shows two snapshots -at different steps of 
degradation- of the lattice in figure 2J The three lat- 
tices (figUJA and fig [5J-B,C) all correspond to a portion 
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FIG. 8: When the angles at junctions are altered by applying 
forces to the edges (as shown in figure [6]), the global orga- 
nization of the pattern is rapidly lost. Correspondingly, the 
slope of the edge rank vs. distance (tg(a)) increases towards 
1. The plot reports this slope as a function of the number of 
iterations of degradations of the original orthogonal Mondrian 
lattice with 2 12 cells. The three letters A, B and C mark the 
values of slope found for the three lattices of which figures [4] 
and [IJ] show a portion. 



of the orthogonal Mondrian lattice with 2 12 cells from 
which the plots of figure [7] are drawn. The plots in figure 
[7J report the statistics of segment histograms (top) and 
edge rank vs. edge distance (bottom) for the three levels 
of degradation of A (no degradation), B (70 iterations 
of degradation) and C (140 iterations). With increas- 
ing degradation, the histogram of segment ranks becomes 
flatter and its center shifts far to the right of the theo- 
retical peak expected for pure Mondrian lattices. The 
slope of edge rank vs. edge distance, that is nearly zero 
for the original lattice A, increases progressively. Figure 
[5] shows that the slope undergoes a continuous transition 
on a relatively large scale (i.e. it is still relatively low also 
for networks, such as the lattice B, whose histogram of 
segment ranks is already flat). This means that once the 
large scale coherence is lost, the length over which co- 
herence persists (inverse of the slope), decreases continu- 
ously. In other words the slope is a quantitative measure 
of the amount of disorder between the two extremes of 
a Mondrian lattice {tg(a) ~ 0) and a foam like pattern 
(tg(a) ~ 1). In the case of this latter, the rank increases 
at almost every junction. 



In general, real world patterns are expected to devi- 
ate in various ways from these simple lattice models. 
In the rest of the paper we explore the distribution of 
ranks and segment statistics in three different examples 
of real world network- like patterns: the pattern of frac- 
ture on the surface of materials, patterns of leaf veins in 
dicotyledon plants, and the pattern of urban streets in 
(unplanned) towns. 



IV. FRACTURE PATTERNS 

Crack patterns often form on the surface of materi- 
als as the result of the shrinking of one material layer 
frustrated by its deposition on a non-shrinking substrate. 
This kind of pattern formation has been extensively ob- 
served and reproduced in controlled settings on a variety 
of materials, including mud, ceramics and coffee grounds. 
The final patterns result from the combination of two 
distinct processes: the nucleation of new fractures and 
the propagation of already existing ones [45l - |47j |. Nucle- 
ation of new fractures usually involves the formation of 
tripartite junctions with equal angles of about 120 de- 
grees |47|- Conversely, junctions formed by propagating 
fractures are the result of either two fractures meeting 
at a point (usually with an orthogonal angle) , or of the 
branching of one growing fracture (in this latter case the 
angles formed at the junction are less predictable). If 
the crack pattern is produced in a non-elastic material, 
there is no reorganization after its formation and the final 
form of the pattern reflects the mechanisms of formation. 
Depending on the characteristics of the material, either 
nucleation of new fractures will be the most frequent pro- 
cess, or propagation of already formed fractures over long 
distances. Nucleation is more frequent in the case of very 
thin layers or inhomogeneous materials; propagation is 
predominant in brittle, homogeneous materials such as 
ceramics. 

Patterns resulting mainly from the propagation of al- 
ready formed fracture lines present a well defined hier- 
archy due to the sequential formation process. One such 
pattern is shown in figure [TJ Figure [S] shows the ranks 
inferred by our algorithm when the two long horizontal 
lines are selected as roots. There is a partial mismatch 
between real and inferred hierarchy. This is in part due 
to the fact that when a cell is cut in two halves by a 
fracture, then the two halves become independent and 
it is no longer possible to establish a temporal relation 
between the new fracture events within a cell and those 




FIG. 9: (Color online) Edge ranks computed for the fracture 
pattern of figure [T] The root edges for rank computation are 
the two horizontal lines that cut the whole figure. The color 
coding is the same as for other figures. The □ and O symbols 
are referred to in the main text. 
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FIG. 10: Crack patterns formed in a layer of paint (top), in thin desiccating clay (middle), in the glaze of ceramics (bottom). 
The graphs on the left plot the histograms of segments with different ranks (as in previous figures); in the case of ceramics, we 
also plot the theoretical distribution for a Mondrian lattice of the same size. The graphs on the central column plot the ranks 
of edges vs. the topological distance of the same edges from the roots. The slope of best fitting straight line is 0.246 for paint, 
0.356 for clay and 0.022 for ceramics. The images on the right are portions of the original patterns, where the fractures have 
been colored according to the rank of the corresponding network edge. The same periodic rainbow colormap as in figure 2] is 
used (color online). The roots for rank computation are all the edges on the left of the picture. 



in the neighboring one. In the case of perfect hierar- 
chical organization (as is the case for the pattern shown 
in figure |9]) one could also gain information by consid- 
ering the organization of angles at both extremities of 
a segment: the fracture line marked with Q terminates 
the fracture marked with □, indicating that this latter is 
actually more recent. Unfortunately, a method that con- 
siders both extremities would lose generality and could 
not resolve a "circular configuration" where a segment 
"A" is terminated by another segment "B" , the segment 
"B" is terminated by "C" and segment "C" is terminated 
by "A". 

Once acknowledged that, depending on the character- 
istics of the material, some fracture patterns are mostly 
dominated by the nucleation of new fractures, and oth- 
ers by the elongation of existing ones, we want to test if 
our algorithm gives different classification results in the 
two cases. To this end, let us consider the crack patterns 
formed in three different materials: paint, desiccating 
clay and ceramics (figure [l"0l from top to bottom). The 
patterns were photographed with a digital camera, con- 
verted to grayscale images, high-pass filtered to remove 
inhomogeneities in the illumination, binarized by simple 



thresholding and cleaned applying a binary morpholog- 
ical majority filter (see [48|, |4t| for a review of common 
image processing techniques). The images of the fracture 
patterns were then skeletonized with a topology preserv- 
ing algorithm based on distance transformation to 
obtain a 8-connected skeleton (a skeleton where two pix- 
els are considered to be connected if they share either a 
face or a corner, opposed to a 4-connected skeleton where 
only pixels that share a face are considered to be con- 
nected). For each skeleton pixel, we counted the number 
of pixels in their 8-neighborhood that also belonged to 
the skeleton, and all the pixels having a number of neigh- 
bors not equal to 2 were marked. All the connected sets 
of marked skeleton pixels were mapped into a network 
node. Whenever there was in the image an unmarked 8- 
connected path between two clusters of pixels identified 
as nodes we introduced an edge between the correspond- 
ing nodes. The orientation of each edge was estimated 
from the coordinates of the two endpoints. In a sub- 
sequent step, edges shorter than a threshold length ( 3 
pixels) were removed and the nodes at the two endpoints 
merged together. 

A rectangular section of the pattern is studied and all 



FIG. 11: (Color online) A Vein pattern of a dicotyledon: Hymenanthera chathamica. A: only the veins with rank and 1 are 
shown in colour (non white lines). B: only veins of rank 0,1 and 2 are colored. C: only veins of rank to 3 are colored. (The 
whole leaf was analyzed, but for clarity only a small portion of the vein network is shown here). 



the edges crossing one side of the rectangle are selected as 
roots for the assignment of ranks. Figure ITU1 reports the 
statistics obtained on the three patterns, together with 
a snapshot of a small region of the original patterns (on 
the right), where the structure is colored according to the 
rank of the corresponding network edge. 

Very different distributions are observed for the cracks 
formed in paint and clay, versus cracks formed in ceram- 
ics: in the former two, the lack of large scale organization 
is reflected into the linear increase of ranks with the topo- 
logical distance from the root and in a flat histogram of 
segment ranks. Indeed, the junctions originating from 
nucleation of new fractures, whose angles are ~ 120 de- 
grees, always determine an increase in the rank of edges 
across the junction. The inverse of the slope of the curve 
fitting the data gives an indication of the length over 
which crack elongation proceeds: about three junctions. 
Conversely, for cracks formed in ceramics, the edge ranks 
are almost completely independent of the topological dis- 
tance of edges from the roots, because of the large scale 
organization typical of these patterns. The slope of the 
rank vs. distance distribution (bottom plot, central col- 
umn) is thus close to zero, but we can still see a deviation 
of the rank histogram toward the left of the theoretical 
normal distribution expected for a pure Mondrian lattice 
(Figure [TU1 bottom left plot). Such small deviation is 
probably due to the fact that the analyzed region is a 
portion of a larger pattern and the cuts introduced by 
the arbitrary frame disconnected some edges from the 
network path that gives them their real rank. This also 
shows that rank distribution is a quite sensitive measure- 
ment. 



V. LEAF VENATION NETWORKS 

Leaf veins in the leaves of flowering plants form char- 
acteristic patterns that can be used by botanists as keys 
for taxonomic identification. However this identification 
is done by eye and does not rely on quantitative measure- 
ments. The pattern is hierarchical, and the diameter of 



a vein roughly reflects its order of appearance during leaf 
morphogenesis, with larger veins being older than smaller 
ones [5l| . Botanists define discrete vein orders looking at 
vein width at the point of branching from its parent vein: 
the large primary vein or midvein is continuous with the 
stem vascular bundles; secondary veins branch from the 
primary vein; tertiary veins are defined by their narrower 
width where they branch from the secondary veins and 
so on. 

We tested our ranking algorithm on the vein pat- 
terns of angiosperm leaves. The leaves were skeletonized 
with 10% sodium hydroxide solution and the pattern was 
scanned with a commercial scanner in transmission mode 
with a resolution superior to 2000 pixels/inch and 256 
gray levels. A network representation was extracted from 
the images in a similar way to what we described for 
cracks. In the ranking procedure, the leaf stem is se- 
lected as root edge. This is in agreement with its special 
role for both transportation and leaf morphogenesis: the 
leaf stem is both the source (through xylem) and the 
sink (through phloem) of all transportation taking place 
in the leaf vein network, as well as the first vein to form 
during leaf morphogenesis [HJ . 

Figure [11] shows a portion of leaf of Hymenanthera 
chathamica with highlighted vein ranks obtained from 
our algorithm: veins with rank equal to 1 or 2 are shown 
in color in A, veins with rank 1 to 3 are highlighted in B 
and veins with rank 1 to 4 in C. The classification does 
not take into account the diameters of the veins, but only 
the orientation of the edges. Nevertheless, the results of 
classification match well the diameter of the veins: from 
the figure, we can see that veins marked with higher ranks 
have in general smaller size (which roughly corresponds 
to say that they have appeared later). In this sense, we 
also recover the vein patterns as derived from the botan- 
ical classification. This is consistent with the hypothesis 
that the older veins are not only larger but also straighter, 
as can be derived from the relation existing between vein 
sizes and junction angles 34]. However, some small veins 
attached to the main vein are given a small rank, while 
one would intuitively ascribe them to a higher rank. This 
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FIG. 12: Top: Segement and edge statistics for Hymenan- 
thera chathamica ( He, top) and Ficus religiosa (Fr, bottom). 
Network size is ~ 12000 nodes for He and ~ 130000 for Fr. 
The slope of the fit for the graphs on the right is 0.019 for He 
and 0.000 For Fr. 

illustrates well the non complete reversibility of the frag- 
mentation process: for each vein we can infer the order 
of appearance with respect to the parent vein, but not 
the exact age. Technically, also the botanical classifica- 
tion -based on diameters- has the same problem, because 
of the existing correlation between branching angles and 
diameter of veins. 

Figure [12] reports the histogram of segment ranks and 
the plot of edge rank vs. distance for both the Hymenan- 
thera chathamica (He) leaf and for a leaf of Ficus reli- 
giosa (Fr). The He network is about one magnitude 
order smaller than the Fr. However, both networks are 
hierarchical: the histogram of segment ranks is peaked 
close to the theoretical value and the fitted slope of edge 
rank vs. distance is close to zero. 

The larger Fr pattern is very close to the perfect hi- 
erarchical pattern (its segment histogram overlaps quite 
well with the one of a perfect Mondrian lattice of the 
same size). Conversely, the histogram of segment ranks 
for He deviates slightly from the theoretical one: it is still 
a Gaussian but shifted to higher values. Unlike the case 
of fractures, this cannot be ascribed to boundary effects 
as the whole leaf with the boundary veins is analysed, 
nor to a simple small number of veins, as the network 
is large and the distribution is already very well defined. 
This shift might be due to the different plasticity of the 
veins during their maturation in the two plant species. 
More measurements should be made to check if this shift 
is due to the difference in species or to the different mat- 
uration of the leaves depending on their sizes. 



VI. URBAN STREET PATTERNS 

The growth of urban street systems is a complex pro- 
cess determined by a series of historical events and involv- 



ing feed-back and regulation from global network pro- 
cesses, such as traffic and transportation. However, in a 
first approximation, for self-organized cities, the general 
form of these patterns can be described in terms of sim- 
ple models where new streets appear over time with no 
reorganisation [13, [EH [Hj] . Within such models, the first 
streets connect the first houses to the country. As the 
urban pattern grows, new streets bifurcate from existing 
ones in the direction of not yet urbanized areas, or to 
join already existing streets. In this respect, the process 
of growth of urban streets shares similarities with the 
growth of fracture patterns (3f| , which supports the fact 
of assigning ranks to streets in a similar way. 

We here test our ranking method on the towns of Cor- 
doba (Spain) and Venice (Italy) and compare the statis- 
tics obtained for the two towns. In both cases we choose 
the perimeter of the town as root for the rank computa- 
tion. (The shores of Venice and the highway ring around 
Cordoba). 

Figure [T3] displays a map of the town of Cordoba, 
where each edge is colored according to its rank. The 
highlighted region in the figure corresponds to the histor- 
ical city centre. Almost all the edges with highest rank 
appear to fall inside this region. Higher ranks often are 
the mark of lack of global organization. A possible out- 
look for future analysis could be testing if higher ranks 
correspond to parts of town that developed in a period of 
more self-organized, organic growth (e.g. periods when 
the central power was weaker). Overall, the histogram 
of segment ranks (reported in figure [T3] top, left) is still 
compatible with that of a hierarchical network, that is, it 
overlaps with the distribution expected for a Mondrian 
lattice of the same size. The distribution of edge ranks 
vs. edge distance from the root edges is nearly flat, with 
a slight slope that can likely be ascribed to the relatively 
small size of the network. 




FIG. 13: (Color online) Cordoba map, edges are colored ac- 
cording to their rank. 
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FIG. 14: Top: Histogram of the number of segments in Cor- 
doba street pattern, Venice street pattern and Venice street 
and channel pattern. Bottom: rank of network edges vs. their 
topological distance from the root edges. The slope of the fit- 
ted line to the bottom plots is 0.094 for Cordoba, 0.239 for 
Venice, and 0.146 when also the channels are considered in 
the network. 



When the same analysis is carried out on the street 
pattern of Venice, a quite different behavior is found: 
segments rank up to much higher values, suggesting an 
absence of large scale organization in this street pattern 
(figure ITU middle column). This does not seem to be 
a pure artifact of our method, but agrees with other el- 
ements of the organization of Venice, where streets are 
not the principal element of organization of the town. 
This is reflected, for instance, in the fact that still to- 
day houses and buildings are not numbered according to 
their position along a street, but following a subdivision 
by districts ("sestieri"). 

Interestingly, when for Venice one considers the trans- 
portation system including both streets and channels, the 
rank distribution gets closer to one of a hierarchical net- 
work (figure [Tfl right). This reveals that the channels 
are the original part of the transportation system, and 
the origin of the town organization. The houses were 
first built along the channels with direct access to them 
ffigure [T5|) , and the streets appeared inside the islands de- 
limited by the channels, as secondary divisions for inland 
house access. Thus removing the channels from the anal- 
ysis is removing the large coherent structure, and thus 
has a direct visible impact on the rank distribution. On 
the contrary, with the channels included, the hierarchi- 
cal structure of successive divisions is recovered. So, the 
rank analysis has the potential of showing some peculiar- 
ities of specific towns, and at the same time of pointing 
to the possible reasons for these peculiarities. 




FIG. 15: Direct access to the channel network from a house 
in Venice 



VII. CONCLUSIONS 

In this article we introduced a new method to ana- 
lyze spatial networks with loops and define quantitative 
measurements on them. The need to characterize quan- 
titatively the form (and the possible morphogenesis) of 
network-like patterns motivated us to introduce a proce- 
dure for assigning ranks to all the edges of the pattern 
and group several edges into segments. Using a local spa- 
tial characteristic, the branching angles, we are able to 
quantify the large scale coherence of the pattern. 

The first measurement we propose, the distributionof 
segment ranks, is a very sensitive indicator of the exis- 
tence of a large scale spatial coherence and of hierarchical 
subdivisions. When this large scale coherence is lost, the 
measurement quantifies the scale over which organiza- 
tion persists, up to the point of purely local organization. 
This scale is directly given by the inverse of the slope of 
edge ranks vs. edge distances from the root. When tested 
on lattice models, these measurements show a dichotomy 
between two types of network, one with large scale co- 
herence and one with purely local organization. 

The method is efficient on different real patterns of 
various origin, from fractures to leaf venation and urban 
streets. The obtained measures show the potential to 
discriminate between patterns with hierarchical history 
of growth and patterns grown out of more local rules. 
For fracture patterns, the two cases of many locally gen- 
erated fractures that reconnect and few fractures prop- 
agating on large distance can be distinguished with no 
ambiguity. For leaf venation pattern, the method clearly 
reveals a hierarchical growth mechanism. In addition, the 
rank assigned to network segments correlates to some ex- 
tent with the temporal order of appearance of the same 
segments, and hence the measure is informative on the 
process of growth itself. The matching between ranks and 
order of appearance however is not perfect, and the infor- 
mation provided by the ranks should be complemented 
from other sources for individual systems. For instance, 
with leaf veins one could consider both ranks and vein 
diameter to obtain a better classification of veins into 
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discrete orders. 

This analysis is also sensitive on town streets, and is 
able to reveal the particular organization of streets in 
Venice, whose structure can be explained as a secondary 
construction from the channels. The town structures are 
otherwise coherent with that of a sequential sub-division 
pattern, indicating some underlying logic in its develop- 
ment. 

We believe that the method can be generally applied 
to different types of patterns, revealing not only their 
structure but also in part the history of their growth. In 
the present paper, we built the classification from local 
angle information. More generally however, any other 
local information can be used, such as -for instance- the 
size of the connecting element, to construct a similar hi- 
erarchy It is thus a new, general and efficient way to 
analyze networks with loops and group patterns of very 
diverse origins into the same structure and development 
classes. 



VIII. APPENDIX 

The distribution of segment ranks in a Mondrian lat- 
tice can be computed more easily if one recognizes that 
the cuts within a cell only affect the distribution of seg- 
ment ranks for that particular cell, but not for the neigh- 
bouring ones. We can identify two kinds of cells ("A" 
and "B" in figure [THJ) , and two possible cuts for each cell 
(either orthogonal or parallel to the side with the low- 
est rank). For any given cell and cut combination, both 
the resulting cells and the new segments added are de- 
termined: a cell of type "A" , whose minimum rank is 
r can give rise either to two cells of type "A" and with 
minimum rank r, or to two cells of type "B" also with 
minimum rank r. In the first case, the cut will introduce 
a new segment with rank r + 1, while in the second it 
will introduce a segment with rank r + 2. A cell of type 



FIG. 16: Schema of the possible cell types in a Mondrian 
lattice and of the results of their divisions. A cell is defined by 
the lowest rank (r) found on one of its sides and by the rank on 
the opposed side to this (the rank of the sides adjacent to that 
with rank r is always r + 1). The cell marked with the letter A 
has rank configuration r, r + 1, r, r + 1 over its perimeter and 
the one marked with B has configuration r, r + 1, r + 2, r + 1. 
For each cell type the cut can be orthogonal to the side with 
rank r or parallel to it. 



"B" and minimum rank r can generate two cells of type 
"B" and rank r; in this case it will add a new segment 
with rank r + 1 as well as a new segment with rank r + 2. 
This latter results from the "cut" operated by the new 
segment with rank r + 1 on the already existing segment 
with rank r + 2 (because we impose that a segment with 
higher rank is always ended when it meets one with lower 
rank) . If the cut is in the other direction, the same "B" 
cell can generate one cell B with minimum rank r and 
one cell "A" with minimum rank r + 1. 

If we assume that all the cuts have the same probabil- 
ity, that all the cells in the lattice are cut once for each 
time step and we take the average realizations, we obtain 
the following equations for the evolution of the pattern: 

A(r,t + l)=A(r,t) + ^B(r-l,t) 

B(r,t+l) = ^B(r,t) + A(r,t) 

Ns(r 7 1 + 1) = Ns(r, t) + B(r -2,t) + ^A(r - 2, t)+ 

+ l -B{r-l,t)+ l -A{r-l,t) 

(4) 

where A(r, t) indicates the number of cells of type "A" 
with minimum rank r at iteration t and similarly for 
B(r,t). Ns{r,t) is the total number of segments with 
rank r at iteration t. 

We iterate the above equations, starting with a single 
cell of type B and minimum rank r — 0. With these 
formulae it is painless to go to very high numbers of cells 
and obtain well defined distributions. The first observa- 
tion in the iterations is that the peaks of the distributions 
of A, B and Ns shifts to the right by one rank unit every 
six iteration steps. These shifts have different phases for 
A, B and Ns. In particular the phase shift for Ns is 
t — 2, that is, the peak position r of Ns increases by one 
unit with pace [ ^"^ j. Figure [171 A plots the position of 
the maximum r as a function of the number of iterations 
for a few iterations of equations |4] around t = 600. 

We also guessed the way the maxima increase in ampli- 
tude with simple relationships: the scale of the numbers 
of cells A, cells B and segments Ns is by contraction 2* 
(the number of cells). However, the maximum value does 
not simply scale in proportion to the number of cells, be- 
cause the distribution gets larger while it moves to the 
right for increasing values of t. In particular, the scaling 
of standard deviation can be easily approximated by a 

We found that if the the 

2 t+i 



simple formula as a = s/t + l. 
distribution is normalized by 



/t+T 



then the maximum 



remains nearly stable through iterations (this is shown 
for a few iterations of equations [4] around t = 600 in 
figure H7J-B). 

When the distribution is rescaled in this way in width, 
position and amplitude, then it falls onto a perfect 
parabola in log coordinates (shown in figure QjJ-C after 
1000 iterations), without any change when increasing the 
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number of iterations, except for attenuating fluctuations. 
So, even if we lack a theoretical demonstration, we know 
that the expression of N in equation |3] is a Guassian and 
the dependencies on the parameters are asymptotically 
exact. Fitting a Gaussian on this distribution after a very 
large number of iterations (t = 1000, that is, a number 
of cells of 2 1000 , -numerical approximation problems ap- 
pear in our code around iteration 1024-), give the fitting 
parameters used in equation [3] of 1.043 and 1.020 respec- 
tively (these numbers could possibly be approximations 
of 7r/3 and respectively). 



Normalized rank (r -2 - 1 (r-2)/6 i)Mr+l) 
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